Quoll Mixed Models Practice
C.S. Stratton
2025-09-29
Intro to Quolls
Does weight of the individual affect max dispersal distance?
GGplot to visualize to start
ggplot(data=data_quoll, aes(Weight, MaxDispersal)) +
geom_point() +
stat_smooth(method="lm", se=FALSE) +
xlab("Weight (g)") +
ylab("Max Dispersal Distance (m)") +
theme_bw() +
theme(legend.title = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())## `geom_smooth()` using formula = 'y ~ x'
Grouped by Sanctuary
ggplot faceted by sanctuary
ggplot(data=data_quoll, aes(Weight, MaxDispersal)) +
geom_point() +
stat_smooth(method="lm", se=FALSE) +
facet_wrap(.~Sanctuary) +
xlab("Weight (g)") +
ylab("Max Dispersal Distance (m)") +
theme_bw() +
theme(legend.title = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())## `geom_smooth()` using formula = 'y ~ x'
model_quoll_sanctuary_1 <- lmer(MaxDispersal ~ Weight + (Weight | Sanctuary), data=data_quoll)## boundary (singular) fit: see help('isSingular')summary(model_quoll_sanctuary_1)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: MaxDispersal ~ Weight + (Weight | Sanctuary)
## Data: data_quoll
##
## REML criterion at convergence: 693
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.2028 -0.6611 -0.2806 0.4212 3.0951
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Sanctuary (Intercept) 3.858e+05 621.1639
## Weight 6.029e-01 0.7765 -1.00
## Residual 4.793e+06 2189.2557
## Number of obs: 39, groups: Sanctuary, 3
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 3420.7662 1548.8921 6.7146 2.209 0.0645 .
## Weight -0.6797 1.4297 2.8904 -0.475 0.6681
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## Weight -0.971
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')Grouped by Sex
ggplot faceted by sex
ggplot(data=data_quoll, aes(Weight, MaxDispersal)) +
geom_point() +
stat_smooth(method="lm", se=FALSE) +
facet_wrap(.~Sex) +
xlab("Weight (g)") +
ylab("Max Dispersal Distance (m)") +
theme_bw() +
theme(legend.title = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank())## `geom_smooth()` using formula = 'y ~ x'
model nested within sex
model_quoll_sex_1 <- lmer(MaxDispersal ~ Weight + (Weight | Sex), data=data_quoll)## boundary (singular) fit: see help('isSingular')summary(model_quoll_sex_1)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: MaxDispersal ~ Weight + (Weight | Sex)
## Data: data_quoll
##
## REML criterion at convergence: 693
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.1655 -0.7117 -0.2735 0.4453 3.0982
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## Sex (Intercept) 2.003e+06 1415.367
## Weight 1.646e+00 1.283 -1.00
## Residual 4.833e+06 2198.412
## Number of obs: 39, groups: Sex, 2
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 3330.51414 1800.11595 0.01682 1.850 0.945
## Weight -0.52157 1.62425 0.01484 -0.321 0.975
##
## Correlation of Fixed Effects:
## (Intr)
## Weight -0.975
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
Will use the sanctuary model going forward
Second Sanctuary Model
Making the model with fixed slopes
model_quoll_sanctuary_2 <- lmer(MaxDispersal ~ Weight + (1 | Sanctuary), data=data_quoll) ## boundary (singular) fit: see help('isSingular')
Making the two flexplots
Random intercepts and random slopes
flexplot::visualize(model_quoll_sanctuary_1, plot="model") + #use all 18 of the measured individuals
ggtitle("Max Dispersal as a function of weight accounting for Sanctuary") +
labs(subtitle = "Formula: MaxDispersal ~ Weight + (Weight | Sanctuary)")## Warning: `aes_string()` was deprecated in ggplot2 3.0.0.
## ℹ Please use tidy evaluation idioms with `aes()`.
## ℹ See also `vignette("ggplot2-in-packages")` for more information.
## ℹ The deprecated feature was likely used in the flexplot package.
## Please report the issue to the authors.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
Random intercepts and fixed slopes
flexplot::visualize(model_quoll_sanctuary_2, plot="model") + #use all 18 of the measured individuals
ggtitle("Max Dispersal as a function of weight accounting for Sanctuary") +
labs(subtitle = "Formula: MaxDispersal ~ Weight + (1 | Sanctuary)")
Display the two models
arm::display(model_quoll_sanctuary_1)## lmer(formula = MaxDispersal ~ Weight + (Weight | Sanctuary),
## data = data_quoll)
## coef.est coef.se
## (Intercept) 3420.77 1548.89
## Weight -0.68 1.43
##
## Error terms:
## Groups Name Std.Dev. Corr
## Sanctuary (Intercept) 621.16
## Weight 0.78 -1.00
## Residual 2189.26
## ---
## number of obs: 39, groups: Sanctuary, 3
## AIC = 705, DIC = 725
## deviance = 709.0arm::display(model_quoll_sanctuary_2)## lmer(formula = MaxDispersal ~ Weight + (1 | Sanctuary), data = data_quoll)
## coef.est coef.se
## (Intercept) 3297.26 1496.11
## Weight -0.56 1.35
##
## Error terms:
## Groups Name Std.Dev.
## Sanctuary (Intercept) 0.00
## Residual 2202.33
## ---
## number of obs: 39, groups: Sanctuary, 3
## AIC = 701.1, DIC = 725
## deviance = 709.0Comparing parameters
compare_quoll_sanctuary <- parameters::compare_parameters(model_quoll_sanctuary_1, model_quoll_sanctuary_2, select = "{estimate}<br>({se})|{p}")
print_html(compare_quoll_sanctuary)| Parameter |
model_quoll_sanctuary_1
|
model_quoll_sanctuary_2
|
||
|---|---|---|---|---|
| Coefficient(SE) | p | Coefficient(SE) | p | |
| (Intercept) | 3420.77 (1548.89) |
0.034 | 3297.26 (1496.11) |
0.034 |
| Weight | -0.68 ( 1.43) |
0.638 | -0.56 ( 1.35) |
0.682 |
| Observations | 39 | 39 | ||
Compare Performance
performance::compare_performance(model_quoll_sanctuary_1, model_quoll_sanctuary_2)## Some of the nested models seem to be identical and probably only vary in
## their random effects.## Random effect variances not available. Returned R2 does not account for random effects.## # Comparison of Model Performance Indices
##
## Name | Model | AIC (weights) | AICc (weights)
## --------------------------------------------------------------------------
## model_quoll_sanctuary_1 | lmerModLmerTest | 721.2 (0.109) | 723.8 (0.056)
## model_quoll_sanctuary_2 | lmerModLmerTest | 717.0 (0.891) | 718.2 (0.944)
##
## Name | BIC (weights) | R2 (cond.) | R2 (marg.) | RMSE
## ----------------------------------------------------------------------------
## model_quoll_sanctuary_1 | 731.2 (0.023) | 0.026 | 0.007 | 2121.649
## model_quoll_sanctuary_2 | 723.7 (0.977) | | 0.004 | 2145.115
##
## Name | Sigma | ICC
## ------------------------------------------
## model_quoll_sanctuary_1 | 2189.256 | 0.019
## model_quoll_sanctuary_2 | 2202.328 |
#### Interaction LM
Comparing perfomance of an
interaction lm just because
model_quoll_sanctuary_lm <- lm(MaxDispersal ~ Weight*Sanctuary, data=data_quoll)
performance::compare_performance(model_quoll_sanctuary_1, model_quoll_sanctuary_2, model_quoll_sanctuary_lm)## Random effect variances not available. Returned R2 does not account for random effects.## # Comparison of Model Performance Indices
##
## Name | Model | AIC (weights) | AICc (weights)
## ---------------------------------------------------------------------------
## model_quoll_sanctuary_1 | lmerModLmerTest | 721.2 (0.083) | 723.8 (0.051)
## model_quoll_sanctuary_2 | lmerModLmerTest | 717.0 (0.683) | 718.2 (0.862)
## model_quoll_sanctuary_lm | lm | 719.2 (0.234) | 722.8 (0.087)
##
## Name | BIC (weights) | RMSE | Sigma | R2 (cond.)
## ---------------------------------------------------------------------------
## model_quoll_sanctuary_1 | 731.2 (0.022) | 2121.649 | 2189.256 | 0.026
## model_quoll_sanctuary_2 | 723.7 (0.951) | 2145.115 | 2202.328 |
## model_quoll_sanctuary_lm | 730.8 (0.027) | 2041.551 | 2219.400 |
##
## Name | R2 (marg.) | ICC | R2 | R2 (adj.)
## -----------------------------------------------------------------
## model_quoll_sanctuary_1 | 0.007 | 0.019 | |
## model_quoll_sanctuary_2 | 0.004 | | |
## model_quoll_sanctuary_lm | | | 0.098 | -0.038